MixVsSingleCompunds

Data loading

suppressPackageStartupMessages({
  library(SummarizedExperiment)  
  library(SEtools)
  library(edgeR)
  library(DT)
  library(pheatmap)
  library(plotly)
  library(dplyr)
  library(sva)
  library(ggrepel)
})
source("Functions/EDC_Functions.R")
source("Functions/CriFormatted.R")
load("Data/AllSEcorrected.RData",verbose=T)

## Loading objects:
##   SEs
##   DEAs

load("Data/HormonalGenes/HormonalGenes.RData",verbose = T)

## Loading objects:
##   Thyroid
##   Androgen
##   Estrogen
##   Corticoid
##   Progesterone
##   PPAR
##   Retinoic

load("Data/DEAfetSingleCompunds.RData", verbose=TRUE)

## Loading objects:
##   res.fetT3
##   res.fetBPA

load("Data/DEAorgSingleCompunds.RData", verbose=TRUE)

## Loading objects:
##   res.orgT3
##   res.orgBPA

MixN vs single compounds Fetal

Fetal MixN DEGs are specific for MixN treated samples

fetal <- SEs$chronic.fetal[,which(SEs$chronic.fetal$EXPO=="CNT" | SEs$chronic.fetal$EXPO=="1X" | SEs$chronic.fetal$EXPO=="1000X" | SEs$chronic.fetal$EXPO=="T3"| SEs$chronic.fetal$EXPO=="BPA0.04X")]

degs.fe <- row.names(DEAs$chronic.fetal)[which((abs(DEAs$chronic.fetal$logFC.EXPO1X)>0.5 | abs(DEAs$chronic.fetal$logFC.EXPO1000X)>0.5) & DEAs$chronic.fetal$FDR<=0.05 & DEAs$chronic.fetal$logCPM>0)]

sehm(fetal, assayName = "corrected", degs.fe,  do.scale = T, show_rownames = F, anno_columns = c("EXPO2","MixtureBatch"), main="MixNDEGs")

Fetal T3 treated samples respond to T3 known targets and the DEGs are specific for T3 samples

sehm(fetal[,which(SEs$chronic.fetal$EXPO=="CNT" | SEs$chronic.fetal$EXPO=="T3")], assayName = "corrected", Thyroid,  do.scale = T, show_rownames = F, anno_columns = c("EXPO2","MixtureBatch"),  main="T3 target genes")

degs.fetT3 <- row.names(res.fetT3)[which(abs(res.fetT3$logFC)> 0.5 & res.fetT3$FDR<=0.05 & res.fetT3$logCPM > 0)]
sehm(fetal, assayName = "corrected", degs.fetT3,  do.scale = T, show_rownames = F, anno_columns = c("EXPO2","MixtureBatch"), main="T3DEGs")

T3 <- res.fetT3
T3$genes <- row.names(T3)

MixN <- DEAs$chronic.fetal
MixN$genes <- row.names(MixN)
MixN$logFC <- (MixN$logFC.EXPO1X + MixN$logFC.EXPO1000X)/2

Comp <- compareResultsFCNew(T3, MixN, FDRth=0.05, FCth=2^0.5, FDRceil=1e-10, logCPMth=0, FCceil=5, 
                              title='Fetal: Comparing T3 to MixN (mean FC)', geneLabel=TRUE, topLab=-30)
Comp$Scatter

#ggsave(Comp$Scatter, filename='Data/SingleCompounds/ScatterFetal_MixN_T3.pdf', width=10, height=10)
#ggsave(Comp$Scatter, filename='Data/SingleCompounds/ScatterFetal_MixN_T3.png', width=10, height=10)

# Check overlap between MixN and Thyroid
table(Comp$AllRes$Status)

## 
## NotSignificant   SignificantA SignificantA&B   SignificantB 
##          12488            374             17             44

## Fisher Test Significant
SigBoth <- dplyr::filter(Comp$Significant, Comp$Significant$Status=='SignificantA&B')

Contingency <- matrix(c(dim(dplyr::filter(SigBoth, log2FCA < 0 & log2FCB < 0 ))[1], 
                          dim(dplyr::filter(SigBoth, log2FCA >= 0 & log2FCB < 0 ))[1],
                           dim(dplyr::filter(SigBoth, log2FCA < 0 & log2FCB >= 0 ))[1], 
                            dim(dplyr::filter(SigBoth, log2FCA >= 0 & log2FCB >= 0 ))[1]),
                        nrow=2, dimnames=list(c('DownA', 'UpA'), c('DownB', 'UpB')))
Contingency

##       DownB UpB
## DownA    11   3
## UpA       1   2

fisher.test(Contingency)

## 
##  Fisher's Exact Test for Count Data
## 
## data:  Contingency
## p-value = 0.1912
## alternative hypothesis: true odds ratio is not equal to 1
## 95 percent confidence interval:
##    0.2518663 468.7960007
## sample estimates:
## odds ratio 
##   6.321431

Fetal BPA treated samples show specific DEGs

degs.fetBPA <- row.names(res.fetBPA)[which(abs(res.fetBPA$logFC)> 0.5 & res.fetBPA$FDR<=0.05 & res.fetBPA$logCPM > 0)]

sehm(fetal, assayName = "corrected", degs.fetBPA,  do.scale = T, show_rownames = F, anno_columns = c("EXPO2","MixtureBatch"), main="BPADEGs")

BPA <- res.fetBPA
BPA$genes <- row.names(BPA)

MixN <- DEAs$chronic.fetal
MixN$genes <- row.names(MixN)
MixN$logFC <- MixN$logFC.EXPO1X

Comp <- compareResultsFCNew(BPA, MixN, FDRth=0.05, FCth=2^0.5, FDRceil=1e-10, logCPMth=0, FCceil=4, 
                              title='Fetal: Comparing BPA to MixN 1X', geneLabel=TRUE, topLab=-30)
Comp$Scatter

#ggsave(Comp$Scatter, filename='Data/SingleCompounds/ScatterFetal_MixN_BPA.pdf', width=10, height=10)
#ggsave(Comp$Scatter, filename='Data/SingleCompounds/ScatterFetal_MixN_BPA.png', width=10, height=10)

MixN vs single compounds Organoids

Organoids MixN DEGs are specific for MixN treated samples

Organoids <- SEs$chronic.org[,which(SEs$chronic.org$EXPO=="CNT"|SEs$chronic.org$EXPO=="1X" | SEs$chronic.org$EXPO=="1000X" | SEs$chronic.org$EXPO=="T3"| SEs$chronic.org$EXPO=="BPA0.04X")]
degs.org <- row.names(DEAs$chronic.org)[which((abs(DEAs$chronic.org$logFC.EXPO1X)>0.5 | abs(DEAs$chronic.org$logFC.EXPO1000X)>0.5) & DEAs$chronic.org$FDR<=0.05 & DEAs$chronic.org$logCPM>0)]

sehm(Organoids, assayName = "corrected", degs.org,  do.scale = T, show_rownames = F, anno_columns = c("EXPO2","MixtureBatch"), main="MixNDEGs")

Organoids T3 treated samples respond to T3 known targets and the DEGs are specific for T3 samples

sehm(Organoids[,which(SEs$chronic.org$EXPO=="CNT" | SEs$chronic.org$EXPO=="T3")], assayName = "corrected", Thyroid,  do.scale = T, show_rownames = F, anno_columns = c("EXPO2","MixtureBatch"), main="T3 target genes")

degs.orgT3 <- row.names(res.orgT3)[which(abs(res.orgT3$logFC)> 0.5 & res.orgT3$FDR<=0.05 & res.orgT3$logCPM > 0)]

sehm(Organoids, assayName = "corrected", degs.orgT3,  do.scale = T, show_rownames = F, anno_columns = c("EXPO2","MixtureBatch"), main="T3DEGs")

T3 <- res.orgT3
T3$genes <- row.names(T3)

MixN <- DEAs$chronic.org
MixN$genes <- row.names(MixN)
MixN$logFC <- (MixN$logFC.EXPO1X + MixN$logFC.EXPO1000X)/2

Comp <- compareResultsFCNew(T3, MixN, FDRth=0.05, FCth=2^0.5, FDRceil=1e-10, logCPMth=0, FCceil=6, 
                              title='Comparing T3 to MixN (mean FC)', geneLabel=TRUE, topLab=-40)
Comp$Scatter

#ggsave(Comp$Scatter, filename='Data/SingleCompounds/ScatterOrg_MixN_T3.pdf', width=10, height=10)
#ggsave(Comp$Scatter, filename='Data/SingleCompounds/ScatterOrg_MixN_T3.png', width=10, height=10)

# Check overlap between MixN and Thyroid
table(Comp$AllRes$Status)

## 
## NotSignificant   SignificantA SignificantA&B   SignificantB 
##          12331           2212             60            286

## Contingency table on shared genes
SigBoth <- dplyr::filter(Comp$Significant, Comp$Significant$Status=='SignificantA&B')
ContingencyII <- matrix(c(dim(dplyr::filter(SigBoth, log2FCA < 0 & log2FCB < 0 ))[1], 
                          dim(dplyr::filter(SigBoth, log2FCA >= 0 & log2FCB < 0 ))[1],
                           dim(dplyr::filter(SigBoth, log2FCA < 0 & log2FCB >= 0 ))[1], 
                            dim(dplyr::filter(SigBoth, log2FCA >= 0 & log2FCB >= 0 ))[1]),
                        nrow=2, dimnames=list(c('DownA', 'UpA'), c('DownB', 'UpB')))
ContingencyII

##       DownB UpB
## DownA    14  36
## UpA      10   0

prop.table(ContingencyII)

##           DownB UpB
## DownA 0.2333333 0.6
## UpA   0.1666667 0.0

fisher.test(ContingencyII)

## 
##  Fisher's Exact Test for Count Data
## 
## data:  ContingencyII
## p-value = 2.601e-05
## alternative hypothesis: true odds ratio is not equal to 1
## 95 percent confidence interval:
##  0.0000000 0.2092738
## sample estimates:
## odds ratio 
##          0

dplyr::filter(SigBoth, log2FCA < 0 & log2FCB > 0 ) %>% dplyr::pull(genes)

##  [1] "SP5"      "COL9A2"   "KRT8"     "DDIT4L"   "IGSF21"   "NEUROG1" 
##  [7] "ID1"      "HIST2H3A" "HIST2H3C" "CDC20B"   "METRN"    "ESM1"    
## [13] "SPATA18"  "RPE65"    "SLC35F2"  "TRIL"     "HIST2H3D" "TRAPPC5" 
## [19] "TMEM121"  "S1PR3"    "SALL3"    "FEZF2"    "PTHLH"    "NDUFS7"  
## [25] "TPH1"     "SOD3"     "ZNF503"   "TGFB2"    "PKMYT1"   "DMRT3"   
## [31] "COL9A3"   "C4orf48"  "CCDC85B"  "APOE"     "COL18A1"  "C19orf24"

Organoids BPA treated samples show specific DEGs

degs.orgBPA <- row.names(res.orgBPA)[which(abs(res.orgBPA$logFC)> 0.5 & res.orgBPA$FDR<=0.05 & res.orgBPA$logCPM > 0)]

sehm(Organoids, assayName = "corrected", degs.orgBPA,  do.scale = T, show_rownames = F, anno_columns = c("EXPO2","MixtureBatch"),  main="BPADEGs")

BPA <- res.orgBPA
BPA$genes <- row.names(BPA)

MixN <- DEAs$chronic.org
MixN$genes <- row.names(MixN)
MixN$logFC <- MixN$logFC.EXPO1X

#Comp <- compareResultsFC(BPA, MixN, FDRth=0.05, FCth=2^0.5, FDRceil=1e-10, FCceil=8,logCPMth = 0, title='Comparing BPA to MixN 1X', geneLabel=TRUE)
#Comp$Scatter

Comp <- compareResultsFCNew(BPA, MixN, FDRth=0.05, FCth=2^0.5, FDRceil=1e-10, logCPMth=0, FCceil=4, 
                              title='Organoids: Comparing BPA to MixN 1X', geneLabel=TRUE, topLab=-30)
Comp$Scatter

#ggsave(Comp$Scatter, filename='Data/SingleCompounds/ScatterOrg_MixN_BPA.pdf', width=10, height=10)
#ggsave(Comp$Scatter, filename='Data/SingleCompounds/ScatterOrg_MixN_BPA.png', width=10, height=10)

# Check overlap between MixN and BPA
table(Comp$AllRes$Status)

## 
## NotSignificant   SignificantA SignificantA&B   SignificantB 
##          13770            821             50            248

## Contingency table on shared genes
SigBoth <- dplyr::filter(Comp$Significant, Comp$Significant$Status=='SignificantA&B')
ContingencyII <- matrix(c(dim(dplyr::filter(SigBoth, log2FCA < 0 & log2FCB < 0 ))[1], 
                          dim(dplyr::filter(SigBoth, log2FCA >= 0 & log2FCB < 0 ))[1],
                           dim(dplyr::filter(SigBoth, log2FCA < 0 & log2FCB >= 0 ))[1], 
                            dim(dplyr::filter(SigBoth, log2FCA >= 0 & log2FCB >= 0 ))[1]),
                        nrow=2, dimnames=list(c('DownA', 'UpA'), c('DownB', 'UpB')))
ContingencyII

##       DownB UpB
## DownA    32   0
## UpA       5  13

prop.table(ContingencyII)

##       DownB  UpB
## DownA  0.64 0.00
## UpA    0.10 0.26

fisher.test(ContingencyII)

## 
##  Fisher's Exact Test for Count Data
## 
## data:  ContingencyII
## p-value = 2.414e-08
## alternative hypothesis: true odds ratio is not equal to 1
## 95 percent confidence interval:
##  12.94172      Inf
## sample estimates:
## odds ratio 
##        Inf

Tentative: overlap

N.B. Overlap analyses with the function below are not equivalent to the scatterplots, because here the selection of DEGs for MixN is done with an ‘or’ approach among the two concentrations for the FC (while above we use MixN 1x FC for BPA and the mean between the two concentrations for T3).

library(overlapper)

## 
## Attaching package: 'overlapper'

## The following object is masked _by_ '.GlobalEnv':
## 
##     overlap.prob

m1 <- overlapper::multintersect(ll=list(fetalDEGs=degs.fe, fetalBPA=degs.fetBPA, fetalT3=degs.fetT3), universe=row.names(DEAs$chronic.fetal), two.tailed = F)

dotplot.multintersect(m1, sizeRange = c(0,15), th=0.05)

## Warning: Removed 1 rows containing missing values (geom_point).

## Warning: Removed 1 rows containing missing values (geom_text).

m2 <- overlapper::multintersect(ll=list(orgDEGs=degs.org, organoidBPA=degs.orgBPA, organoidT3=degs.orgT3), universe=row.names(DEAs$chronic.org),two.tailed = F)

dotplot.multintersect(m2, sizeRange = c(0,15), th=0.05)

## Warning: Removed 1 rows containing missing values (geom_point).

## Warning: Removed 1 rows containing missing values (geom_text).

Data praparation

For details on data filtering, normalization, batch correction and differential expression analysis, see here and here

---

Index

MixN_Acute

MixN_ChronicFetal

FunctionalFetalChronicMixN

MixN_ChronicOrganoids

FunctionalOrganoidChronicMixN

MixNVsSingleCompounds

MasterRegulatorAnalysis

HormonalPathways

Immunofluorescence

MixN_clustering

FunctionalClustering

Authors

Nicolò Caporale: nicolo.caporale@ieo.it, nicolo.caporale@unimi.it

Cristina Cheroni: cristina.cheroni@ieo.it

Pierre-Luc Germain: pierre-luc.germain@ieo.it

Giuseppe Testa: giuseppe.testa@ieo.it

Lab: http://www.testalab.eu/

‘Date: April 09, 2021’

sessionInfo()

## R version 4.0.3 (2020-10-10)
## Platform: x86_64-apple-darwin17.0 (64-bit)
## Running under: macOS Big Sur 10.16
## 
## Matrix products: default
## BLAS:   /Library/Frameworks/R.framework/Versions/4.0/Resources/lib/libRblas.dylib
## LAPACK: /Library/Frameworks/R.framework/Versions/4.0/Resources/lib/libRlapack.dylib
## 
## locale:
## [1] en_US.UTF-8/en_US.UTF-8/en_US.UTF-8/C/en_US.UTF-8/en_US.UTF-8
## 
## attached base packages:
## [1] parallel  stats4    stats     graphics  grDevices utils     datasets 
## [8] methods   base     
## 
## other attached packages:
##  [1] overlapper_0.99.1           ggrepel_0.9.1              
##  [3] sva_3.38.0                  BiocParallel_1.24.1        
##  [5] genefilter_1.72.1           mgcv_1.8-34                
##  [7] nlme_3.1-152                dplyr_1.0.5                
##  [9] plotly_4.9.3                ggplot2_3.3.3              
## [11] pheatmap_1.0.12             DT_0.17                    
## [13] edgeR_3.32.1                limma_3.46.0               
## [15] SEtools_1.4.0               SummarizedExperiment_1.20.0
## [17] Biobase_2.50.0              GenomicRanges_1.42.0       
## [19] GenomeInfoDb_1.26.2         IRanges_2.24.1             
## [21] S4Vectors_0.28.1            BiocGenerics_0.36.0        
## [23] MatrixGenerics_1.2.1        matrixStats_0.58.0         
## 
## loaded via a namespace (and not attached):
##   [1] circlize_0.4.12        plyr_1.8.6             lazyeval_0.2.2        
##   [4] shinydashboard_0.7.1   splines_4.0.3          digest_0.6.27         
##   [7] foreach_1.5.1          htmltools_0.5.1.1      fansi_0.4.2           
##  [10] magrittr_2.0.1         memoise_2.0.0          cluster_2.1.1         
##  [13] openxlsx_4.2.3         ComplexHeatmap_2.6.2   annotate_1.68.0       
##  [16] colorspace_2.0-0       blob_1.2.1             xfun_0.21             
##  [19] crayon_1.4.1           RCurl_1.98-1.2         jsonlite_1.7.2        
##  [22] survival_3.2-7         iterators_1.0.13       glue_1.4.2            
##  [25] registry_0.5-1         gtable_0.3.0           zlibbioc_1.36.0       
##  [28] XVector_0.30.0         UpSetR_1.4.0           GetoptLong_1.0.5      
##  [31] DelayedArray_0.16.1    V8_3.4.0               shape_1.4.5           
##  [34] scales_1.1.1           futile.options_1.0.1   DBI_1.1.1             
##  [37] randomcoloR_1.1.0.1    Rcpp_1.0.6             viridisLite_0.3.0     
##  [40] xtable_1.8-4           clue_0.3-58            bit_4.0.4             
##  [43] htmlwidgets_1.5.3      httr_1.4.2             RColorBrewer_1.1-2    
##  [46] ellipsis_0.3.1         pkgconfig_2.0.3        XML_3.99-0.5          
##  [49] farver_2.1.0           sass_0.3.1             locfit_1.5-9.4        
##  [52] utf8_1.2.1             reshape2_1.4.4         tidyselect_1.1.0      
##  [55] labeling_0.4.2         rlang_0.4.10           later_1.1.0.1         
##  [58] AnnotationDbi_1.52.0   munsell_0.5.0          tools_4.0.3           
##  [61] cachem_1.0.4           generics_0.1.0         RSQLite_2.2.3         
##  [64] evaluate_0.14          stringr_1.4.0          fastmap_1.1.0         
##  [67] yaml_2.2.1             knitr_1.31             bit64_4.0.5           
##  [70] shinycssloaders_1.0.0  zip_2.1.1              purrr_0.3.4           
##  [73] mime_0.10              formatR_1.7            compiler_4.0.3        
##  [76] curl_4.3               png_0.1-7              tibble_3.1.0          
##  [79] bslib_0.2.4            stringi_1.5.3          futile.logger_1.4.3   
##  [82] highr_0.8              lattice_0.20-41        Matrix_1.3-2          
##  [85] ggsci_2.9              vctrs_0.3.7            pillar_1.5.1          
##  [88] lifecycle_1.0.0        jquerylib_0.1.3        GlobalOptions_0.1.2   
##  [91] cowplot_1.1.1          data.table_1.14.0      bitops_1.0-6          
##  [94] seriation_1.2-9        httpuv_1.5.5           R6_2.5.0              
##  [97] promises_1.2.0.1       TSP_1.1-10             gridExtra_2.3         
## [100] codetools_0.2-18       lambda.r_1.2.4         assertthat_0.2.1      
## [103] rjson_0.2.20           withr_2.4.1            GenomeInfoDbData_1.2.4
## [106] VennDiagram_1.6.20     grid_4.0.3             tidyr_1.1.3           
## [109] rmarkdown_2.7          Cairo_1.5-12.2         Rtsne_0.15            
## [112] shiny_1.6.0